Adding operator support to time evolution problem #606

Gavin-Rockwood · 2025-12-01T05:49:16Z

Checklist

Thank you for contributing to QuantumToolbox.jl! Please make sure you have finished the following tasks before opening the PR.

Please read Contributing to Quantum Toolbox in Julia.
Any code changes were done in a way that does not break public API.
Appropriate tests were added and tested locally by running: make test.
Any code changes should be julia formatted by running: make format.
All documents (in docs/ folder) related to code changes were updated and able to build locally by running: make docs.
(If necessary) the CHANGELOG.md should be updated (regarding to the code changes) and built by running: make changelog.

Request for a review after you have completed all the tasks. If you have not finished them all, you can also open a Draft Pull Request to let the others know this on-going work.

Description

Updated TimeEvolutionProblem be parametrically typed using the type of the starting state as well as expanded sesolve, mesolve and mcsolve to allow them to take in operators, superoperators and operators respectively (along with other changes to facilitate this). ssesolve and smesolve are still in progress.

Related issues or PRs

Working on #521

…nsions for examples.

…ions for usage.

…onProblem

…ense.

.gitignore

src/time_evolution/time_evolution.jl

src/time_evolution/mcsolve.jl

src/time_evolution/mesolve.jl

src/time_evolution/mcsolve.jl

…onProblem

…ussed in the PR comments.

src/time_evolution/mcsolve.jl

ytdHuang · 2025-12-04T14:46:00Z

src/time_evolution/mcsolve.jl

    output_func::Union{Tuple,Nothing} = nothing,
    kwargs...,
-) where {TJC<:LindbladJumpCallbackType}
+) where {TJC<:LindbladJumpCallbackType,ST<:Union{Ket,Operator}}


Suggested change

) where {TJC<:LindbladJumpCallbackType,ST<:Union{Ket,Operator}}

) where {ST<:Union{Ket,Operator},TJC<:LindbladJumpCallbackType}

I actually put this back to just ST<:Union{Ket}. mcsolve returns the average density matrices of all trajectories. Extra work would need to be done to expand _average_traj_states to take a matrix of operators, turn each column into kets and then averaging the density matrix for each of those kets for all trajectories. I don't think this would be hard to do but I'm not sure how useful it would be?

ytdHuang · 2025-12-04T14:47:02Z

src/time_evolution/mcsolve.jl

    ensemble_prob = TimeEvolutionProblem(
        EnsembleProblem(prob_mc.prob, prob_func = _prob_func, output_func = _output_func[1], safetycopy = false),
        prob_mc.times,
+        ST(),


Suggested change

ST(),

ψ0.type,

This was replaced with prob_mc.states_type.

ytdHuang · 2025-12-04T14:47:19Z

src/time_evolution/mcsolve.jl

    normalize_states::Union{Val,Bool} = Val(true),
    kwargs...,
-) where {TJC<:LindbladJumpCallbackType}
+) where {TJC<:LindbladJumpCallbackType} where {ST<:Union{Ket,Operator}}


Suggested change

) where {TJC<:LindbladJumpCallbackType} where {ST<:Union{Ket,Operator}}

) where {ST<:Union{Ket,Operator},TJC<:LindbladJumpCallbackType}

ytdHuang · 2025-12-04T14:50:55Z

src/time_evolution/mesolve.jl

    ens_prob = TimeEvolutionProblem(
        EnsembleProblem(prob.prob, prob_func = _prob_func, output_func = _output_func[1], safetycopy = false),
        prob.times,
+        StateOpType(),


Suggested change

StateOpType(),

prob.state_type,

ytdHuang · 2025-12-04T14:51:31Z

src/time_evolution/sesolve.jl

    prob = ODEProblem{getVal(inplace),FullSpecialize}(U, ψ0, tspan, params; kwargs4...)

-    return TimeEvolutionProblem(prob, tlist, H_evo.dimensions)
+    return TimeEvolutionProblem(prob, tlist, ST(), H_evo.dimensions)


Suggested change

return TimeEvolutionProblem(prob, tlist, ST(), H_evo.dimensions)

return TimeEvolutionProblem(prob, tlist, ψ0.type, H_evo.dimensions)

This is needed because psi0 gets overwritten in line 82 and is no longer a quantum object in line 93. I can leave like this or I can put a line states_type = ѱ0.data somewhere before 82.

Edit: actually that is exactly what I did for mesolve so I will do the same for sesolve and elsewhere.

ytdHuang · 2025-12-04T14:53:02Z

src/time_evolution/smesolve.jl

    )

-    return TimeEvolutionProblem(prob, tlist, dims, (isoperket = Val(isoperket(ψ0)),))
+    return TimeEvolutionProblem(prob, tlist,StateOpType(), dims, ())


Suggested change

return TimeEvolutionProblem(prob, tlist,StateOpType(), dims, ())

return TimeEvolutionProblem(prob, tlist,ψ0.type, dims)

Similar to mesolve, I introduced a states_type variable that is defined during the definition of $\rho 0$.

ytdHuang · 2025-12-04T14:54:00Z

src/time_evolution/smesolve.jl

    ensemble_prob = TimeEvolutionProblem(
        EnsembleProblem(prob_sme, prob_func = _prob_func, output_func = _output_func[1], safetycopy = true),
        prob_sme.times,
+        StateOpType(),


Suggested change

StateOpType(),

ψ0.type,

this was replaced with prob_sme.states_type.

ytdHuang · 2025-12-04T14:54:43Z

src/time_evolution/ssesolve.jl

    )

-    return TimeEvolutionProblem(prob, tlist, dims)
+    return TimeEvolutionProblem(prob, tlist, ST(), dims)


Suggested change

return TimeEvolutionProblem(prob, tlist, ST(), dims)

return TimeEvolutionProblem(prob, tlist, ψ0.type, dims)

Similar to sesolve, this got updated to a states_type variable defined near the top of the function.

ytdHuang · 2025-12-04T14:56:00Z

src/time_evolution/ssesolve.jl

    ensemble_prob = TimeEvolutionProblem(
        EnsembleProblem(prob_sme, prob_func = _prob_func, output_func = _output_func[1], safetycopy = true),
        prob_sme.times,
+        ST(),


Suggested change

ST(),

ψ0.type,

This became prob_sme.times. (actually prob_sse, I changed it to be consistent). Though, as of right now, the only type ssesolve will take is just Ket. The there is some issue with how the problem is defined that makes it so it cannot take in matrices. Same also didn't update smesolve for the same reason.

ytdHuang · 2025-12-06T03:45:27Z

src/time_evolution/mesolve.jl

    T = Base.promote_eltype(L_evo, ψ0)
-    ρ0 = if isoperket(ψ0) # Convert it to dense vector with complex element type
-        to_dense(_complex_float_type(T), copy(ψ0.data))
-    else
-        to_dense(_complex_float_type(T), mat2vec(ket2dm(ψ0).data))
+    # ρ0 = if isoperket(ψ0) # Convert it to dense vector with complex element type
+    #     to_dense(_complex_float_type(T), copy(ψ0.data))
+    # else
+    #     to_dense(_complex_float_type(T), mat2vec(ket2dm(ψ0).data))
+    # end
+    if isoper(ψ0)
+        ρ0 = to_dense(_complex_float_type(T), mat2vec(ψ0.data))
+        state_type = Operator()
+    elseif isoperket(ψ0)
+        ρ0 = to_dense(_complex_float_type(T), copy(ψ0.data))
+        state_type = OperatorKet()
+    elseif isket(ψ0)
+        ρ0 = to_dense(_complex_float_type(T), mat2vec(ket2dm(ψ0).data))
+        state_type = Operator()
+    elseif issuper(ψ0)
+        ρ0 = to_dense(_complex_float_type(T), copy(ψ0.data))
+        state_type = SuperOperator()
    end


Suggested change

T = Base.promote_eltype(L_evo, ψ0)

ρ0 = if isoperket(ψ0) # Convert it to dense vector with complex element type

to_dense(_complex_float_type(T), copy(ψ0.data))

else

to_dense(_complex_float_type(T), mat2vec(ket2dm(ψ0).data))

# ρ0 = if isoperket(ψ0) # Convert it to dense vector with complex element type

# to_dense(_complex_float_type(T), copy(ψ0.data))

# else

# to_dense(_complex_float_type(T), mat2vec(ket2dm(ψ0).data))

# end

if isoper(ψ0)

ρ0 = to_dense(_complex_float_type(T), mat2vec(ψ0.data))

state_type = Operator()

elseif isoperket(ψ0)

ρ0 = to_dense(_complex_float_type(T), copy(ψ0.data))

state_type = OperatorKet()

elseif isket(ψ0)

ρ0 = to_dense(_complex_float_type(T), mat2vec(ket2dm(ψ0).data))

state_type = Operator()

elseif issuper(ψ0)

ρ0 = to_dense(_complex_float_type(T), copy(ψ0.data))

state_type = SuperOperator()

end

# Convert to dense vector with complex element type

T = _complex_float_type(Base.promote_eltype(L_evo, ψ0))

if isoperket(ψ0) || issuper(ψ0)

ρ0 = to_dense(T, copy(ψ0.data))

states_type = ψ0.type

else

ρ0 = to_dense(T, mat2vec(ket2dm(ψ0).data))

states_type = Operator()

end

ket2dm will do nothing if the input is density matrix (Operator), so it doesn't hurt performance, and we can also do ket2dm for Operator input.

Also, maybe also use states_type instead of state_type to avoid confusion, since we call it states_type in TimeEvolutionProblem

ytdHuang · 2025-12-06T03:53:49Z

src/time_evolution/mesolve.jl

+        if isoper(state)
+            to_dense(_complex_float_type(T), mat2vec(state.data))
+        elseif isoperket(state)
+            to_dense(_complex_float_type(T), copy(state.data))
+        elseif isket(state)
+            to_dense(_complex_float_type(T), mat2vec(ket2dm(state).data))
+        elseif issuper(state)
+            to_dense(_complex_float_type(T), copy(state.data))
        end


same here

Suggested change

if isoper(state)

to_dense(_complex_float_type(T), mat2vec(state.data))

elseif isoperket(state)

to_dense(_complex_float_type(T), copy(state.data))

elseif isket(state)

to_dense(_complex_float_type(T), mat2vec(ket2dm(state).data))

elseif issuper(state)

to_dense(_complex_float_type(T), copy(state.data))

end

if isoperket(state) || issuper(state)

to_dense(T, copy(state.data))

else

to_dense(T, mat2vec(ket2dm(state).data))

end

ytdHuang · 2025-12-06T03:55:49Z

src/time_evolution/mesolve.jl

-    if getVal(isoperket)
-        ρt = map(ϕ -> QuantumObject(ϕ, type = OperatorKet(), dims = dimensions), sol.u)
+function _gen_mesolve_solution(sol, prob::TimeEvolutionProblem{ST}) where {ST<:Union{Operator,OperatorKet,SuperOperator}}
+    if prob.states_type == Operator


I think here should be

Suggested change

if prob.states_type == Operator

if prob.states_type isa Operator

cause prob.states_type should be Operator() ?

Gavin-Rockwood added 5 commits November 30, 2025 22:25

made most of the changes. Look at testing.ipynb ##Testing Solver Exte…

571abc9

…nsions for examples.

Added operator support. Look at testing.ipynb ##Testing Solver Extens…

9665194

…ions for usage.

Merge branch 'qutip:main' into Adding-Operator-Support-to-TimeEvoluti…

0e6036a

…onProblem

updated testing.ipynb

01baa16

added to_dense call to sesolve_map to ensure all initial states are d…

e575f7d

…ense.

ytdHuang reviewed Dec 1, 2025

View reviewed changes

Gavin-Rockwood added 4 commits December 2, 2025 18:40

fixed .gitignore and removed testing.ipynb

71842a7

Merge branch 'qutip:main' into Adding-Operator-Support-to-TimeEvoluti…

41e8216

…onProblem

Updated the state type variable to ST and made the other changes disc…

5a354fe

…ussed in the PR comments.

fixed a bug with smesolve!

d8d26a3

ytdHuang reviewed Dec 4, 2025

View reviewed changes

ytdHuang reviewed Dec 6, 2025

View reviewed changes

	) where {TJC<:LindbladJumpCallbackType,ST<:Union{Ket,Operator}}
	) where {ST<:Union{Ket,Operator},TJC<:LindbladJumpCallbackType}

	) where {TJC<:LindbladJumpCallbackType} where {ST<:Union{Ket,Operator}}
	) where {ST<:Union{Ket,Operator},TJC<:LindbladJumpCallbackType}

	return TimeEvolutionProblem(prob, tlist, ST(), H_evo.dimensions)
	return TimeEvolutionProblem(prob, tlist, ψ0.type, H_evo.dimensions)

	if prob.states_type == Operator
	if prob.states_type isa Operator

Adding operator support to time evolution problem #606

Are you sure you want to change the base?

Adding operator support to time evolution problem #606

Conversation

Gavin-Rockwood commented Dec 1, 2025 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Checklist

Description

Related issues or PRs

Uh oh!

Uh oh!

Uh oh!

Uh oh!

Uh oh!

Uh oh!

Uh oh!

Uh oh!

Uh oh!

Uh oh!

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Gavin-Rockwood Dec 4, 2025 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Gavin-Rockwood Dec 4, 2025 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Reviewers

Assignees

Labels

Projects

Milestone

Development

Uh oh!

2 participants

Gavin-Rockwood commented Dec 1, 2025 •

edited

Loading

Gavin-Rockwood Dec 4, 2025 •

edited

Loading

Gavin-Rockwood Dec 4, 2025 •

edited

Loading